import numpy as np
import pandas as pd
import geatpy as ea
"""
该程序用于第三问的优化求解：
主要算法方法：NSGA-II（Non-dominated Sorting Genetic Algorithm II）
该算法主要用于多目标规划
改进之处：
1 融入初始值、引入全局存档加快模型收敛
2 引入了混合编码（'BG' 编码-解决离散的，'RI' 编码-解决连续的）
"""
address = r'data/q3.csv'
data = pd.read_csv(address)
data_40_rows = data.head(40)

cost_pre = data_40_rows['批发价格(元/千克)'].tolist() # 预测成本
sale_pre = data_40_rows['平均销量'].tolist() # 预测销量
for i in range(40):
    if sale_pre[i] <= 2.5:
        sale_pre[i] = 2.5
mean_sale = data_40_rows['平均单价']
alph = [0.1]* 40 # 损耗率
discount = [0.8] * 40 # 打折
ini_data = np.array([1]*30 + [0]*10 + [0.7]*40 + sale_pre) # 初始值
need_list = data_40_rows['Q'].tolist()
weight = [0.5] * 40

#销量修正函数，与遗传无关
def Modify(sale_price):
    sale = []
    for i in range(40):
        ss = ( 1 - (sale_price[i] - mean_sale[i])/(mean_sale[i]) * weight[i] ) + sale_pre[i]
        sale.append(ss)
    return sale

#目标函数
def Obj_fuction(X):
    Objv= []#目标值数组
    Cv = []#约束矩阵数组
    #对种群每一个遍历
    for i in range(len(X)):
        x = X[i]
        # 目标一：最大化利润
        decide = x[:40] #01变量
        profit = x[40:80] #利润率
        nums = x[80:120] #补货量
        # 计算成本
        cost = []
        for i in range(len(cost_pre)):
            cost.append(decide[i] * cost_pre[i])
        # 计算售价
        sale_price = []
        for i in range(len(cost_pre)):
            s = cost[i] * (1 + profit[i])
            sale_price.append(s)
        # 计算销量
        sale_modify = Modify(sale_price)
        # 计算总销售额
        w1 = []
        for i in range(len(sale_modify)):
            bad = (sale_modify[i] * alph[i]) * ( discount[i] * sale_price[i] )
            good = (sale_modify[i] * (1 - alph[i])) *(sale_price[i])
            w1.append(bad + good)
        sum_sale = sum(w1)
        # 计算总成本
        w2 = []
        for i in range(len(cost)):
            w2.append(nums[i] * cost[i])
        sum_cost = sum(w2)
        # 惩罚函数
        penalty = 0
        for i in range(40):
            if nums[i] < 2.5:
                penalty -= 200
        f1 = sum_sale - sum_cost + penalty
        # 目标二：最大化需求
        satif_need = []
        for i in range(len(decide)):
            satif_need.append(decide[i] * need_list[i])
        f2 = sum(satif_need) / sum(need_list)

        # 确定评价函数值
        x_Objv = np.hstack([f1,f2]) #目标函数矩阵
        # 定义约束
        list_cv = [sum(decide) - 33, 27 - sum(decide)] #违反约束矩阵：约束条件<=0
        # for i in range(40):
        #   list_cv.append( 2.5- nums[i] )
        x_Cv = np.array(list_cv)
        Objv.append(x_Objv)
        Cv.append(x_Cv)
    Objv = np.array(Objv)
    Cv = np.array(Cv)
    return Objv, Cv
problem = ea.Problem(
    name='NSGAII_for_q3',
    M = 2,#目标维数
    maxormins = [-1,-1],#1是最小化，-1是最大化
    Dim = 120,#决策变量维数
    varTypes=[1] * 40 + [0] * 80,#决策变量列表类型，0：实数，1：整数
    lb=[0] * 40 + [0] * 40 + data_40_rows['最小销量'].tolist(),#下界
    ub=[1] * 40 + [1.2] * 40 + data_40_rows['最大销量'].tolist(),#上界
    evalVars=Obj_fuction
)
# 混合编码
Encodings=['BG','RI']#BG是二进制编码，RI是实数编码
Field1 = ea.crtfld( Encodings[0], problem.varTypes[:40],ranges=np.array( [problem.lb[:40],problem.ub[:40]] ) )#定义前40个的编码解码
Field2 = ea.crtfld( Encodings[1], problem.varTypes[40:],ranges=np.array( [problem.lb[40:],problem.ub[40:]] ) )#定义40以后的编码解码
Fields = [Field1, Field2]
population = ea.PsyPopulation(Encodings=Encodings, Fields=Fields, NIND=100)#NIND是初始种群数

algorithm = ea.moea_psy_NSGA2_archive_templet(
    problem,
    population,
    MAXGEN=300,#最大迭代次数
    logTras=20,#每20次追踪一次
    prophetPop=ini_data,#预设的先知种群
    maxTrappedCount=20,#最大陷入计数避免局部最优
)
res = ea.optimize(algorithm, seed=1, verbose=True , drawing=1,outputMsg=True, drawLog = True )
print(f" 最优解是：{res['Vars'][0]}")
print(f" 最优解值是：{res['ObjV'][0]}")
print(f"hv 结果为：{res['hv']}")
print(f"spacing 结果为：{res['spacing']}")
'''
输出参数:
一个保存着结果的字典。内容为：
    {'success': True or False,  # 表示算法是否成功求解。
    'stopMsg': xxx,  # 存储着算法停止原因的字符串。
    'optPop': xxx,  # 存储着算法求解结果的种群对象。如果无可行解，则optPop.sizes=0。optPop.Phen为决策变量矩阵，optPop.ObjV为目标函数值矩阵。
    'lastPop': xxx,  # 算法进化结束后的最后一代种群对象。
    'Vars': xxx,  # 等于optPop.Phen，此处即最优解。若无可行解，则Vars=None。
    'ObjV': xxx,  # 等于optPop.ObjV，此处即最优解对应的目标函数值。若无可行解，ObjV=None。
    'CV': xxx,  # 等于optPop.CV，此处即最优解对应的违反约束程度矩阵。若无可行解，CV=None。
    'startTime': xxx,  # 程序执行开始时间。
    'endTime': xxx,  # 程序执行结束时间。
    'executeTime': xxx,  # 算法所用时间。
    'nfev': xxx,  # 算法评价次数
    'gd': xxx,  # (多目标优化且给定了理论最优解时才有) GD指标值。
    'igd': xxx,  # (多目标优化且给定了理论最优解时才有) IGD指标值。
    'hv': xxx,  # (多目标优化才有) HV指标值。
    'spacing': xxx}  # (多目标优化才有) Spacing指标值。}
'''